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We formulate a theory of transport in graphene bilayers in the weak momentum scattering regime 
in such a way as to take into account contributions to the electrical conductivity to leading and 
next-to-leading order in the scattering potential. The response of bilayers to an electric field cannot 
be regarded as a sum of terms due to individual layers. Rather, interlayer tunneling and coherence 
between positive- and negative-energy states give the main contributions to the conductivity. At low 
energies, the dominant effect of scattering on transport comes from scattering within each energy 
band, yet a simple picture encapsulating the role of collisions in a set of scattering times is not 
applicable. Coherence between positive- and negative-energy states gives, as in monolayers, a term 
in the conductivity which depends on the order of limits. The application of an external gate, 
which introduces a gap between positive- and negative-energy states, does not affect transport. 
Nevertheless the solution to the kinetic equation in the presence of such a gate is very revealing for 
transport in both bilayers and monolayers. 



I. INTRODUCTION 

The manufacture of carbon monolayers in the 
laboratory has generated unprecedented excitement. 
The unusual and sometimes baffling properties of 
graphene have exposed the community to novel sci- 
ence and supplied ideas for technological innovation. 1 
This achievement was swiftly followed by the reliable 
manufacture of graphene bilayers and multilayers, 
which have become independent research areas in 
their own right. The electronic properties of bilayer 
graphene 2 - have been the subject of recent reviews 3 -^ 
but the topic has been studied much less than single 
graphene layers. Large mobilities recently observed^ 
are among many factors justifying experimental and 
theoretical interest in bilayers. Also notable in these 
systems is a Berry phase of 2-7T accompanying an 
unconventional quantum Hall effect^ and predictions 
of Andreev reflection 7 - and superfluidity^ among other 
effects^ilLi^il^^I Studies of graphene bilay- 
ers have focused on transport f 18 ' 19 ' 20 ' 21 ' 22 ' 23 ' 24 ' 25 ' 26 ' 27 



electron-electron 
and electronic 



compressibility^ 8 - impurities ) 29 ! 30 ' 31 
interactions ) 32 ' 33 ! 3 ^^ and band 
structure^ 37 ! 38 ! 39 ! 40 - 41 

Graphene bilayers are interesting from a physics point 
of view because they are not merely the sum of two lay- 
ers. Tunneling between layers is characterized by the 
large parameter t± which is comparable to the Fermi 
energy at carrier densities n well beyond 10 12 cm~ 2 . 
The energy spectrum of bilayers consists of four bands, 
two with positive energy and two with negative energy. 
The lower positive-energy band and the higher negative- 
energy band touch at k — 0. Although at low densities 
only one of these latter two bands is occupied, depending 
on whether the sample is doped with electrons or holes, 
a pseudospin cannot be defined in the same manner as 
in single graphene layers. This is due to the presence 
of the large tunneling parameter, which indicates a non- 
trivial interplay of carriers from both layers, particularly 



important in steady-state processes. Furthermore, a gate 
potential opens a gap in the energy spectrum^ and in- 
dependent control of this gap and of the carrier density 
can be achieved by using separate top and back gates42 

This article presents a thorough investigation of trans- 
port in graphene bilayers, constructing a compact and 
straightforward framework for analyzing the structure of 
the steady-state density matrix in an electric field. The 
formalism takes the density operator and quantum Li- 
ouville equation as its starting point, treating all terms 
in the Hamiltonian on the same footing and using real- 
istic scattering potentials. Calculation of the electrical 
current reveals the complex physics underlying bilayer 
transport. The conductivity contains a term which is a 
function of the carrier density n and is inversely propor- 
tional to the impurity density rii, similar to the usual 
conductivity of metals and semiconductors. This term 
also depends on t± , which can be thought of as indicat- 
ing coherence between layers, while its dependence on n 
is different for long-range and short range impurities and 
it is in general not oc n. The concept of a characteristic 
momentum scattering time is useful as an order of mag- 
nitude, but it is not possible to assign scattering times 
to carriers in different bands. Furthermore, an electric 
field couples positive and negative energy states, result- 
ing in a term in the conductivity similar to that in single 
layers but of different magnitude because of the differ- 
ent winding number associated with energy dispersion 
in bilayers. ft must also be borne in mind that positive- 
and negative-energy states involve carriers from both lay- 
ers. In graphene monolayers the conductivity due to the 
coupling of positive- and negative-energy states is renor- 
malized by scattering to twice its original value. In bi- 
layers the effect of scattering on this term is considerably 
smaller than in monolayers. This term also depends on 
the order of limits, and work to date has not been able 
to extract its definite value. Fortunately, we will show 
that biased bilayer graphene ought to provide an answer. 
The external gate makes a contribution to the conduc- 
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tivity through an off-diagonal (Hall) term which would 
appear to exist without a magnetic field. This term also 
depends on the order of limits, but it is the sole contribu- 
tion to the off-diagonal conductivity. Since crystal sym- 
metry and Onsager relations imply that an off-diagonal 
conductivity cannot exist without a magnetic field this 
must indicate the correct and unambiguous order of lim- 
its. If any doubts persist experiment can surely resolve 
this issue, given that there is only one potential contri- 
bution to the off-diagonal conductivity. We note that 
Boltzmann transport theory has been formulated thor- 
oughly for unbiased bilayer graphene^ In this work we 
wish to consider both the ordinary ( "Boltzmann" ) contri- 
bution to the conductivity oc n/rii and the contribution 
independent of n and rn on the same footing. We also 
aim to determine unambiguously the role of the gate and 
the lessons to be learned from it. 

We focus on extrinsic graphene in the weak momen- 
tum scattering regime, ep 3> H/t, where ep is the Fermi 
energy and r cx n^ 1 is a characteristic momentum scat- 
tering time (for strong momentum scattering ep <C S/t.) 
Extrinsic graphene refers to the case n > and may be 
in either the weak or the strong momentum scattering 
regime. Intrinsic graphene refers to the case n = where 
the Fermi energy lies at the point where the bands touch, 
and is by definition in the strong momentum scattering 
regime. Enormous strides in sample quality make trans- 
port in the weak momentum scattering regime a timely 
undertaking^*^ We assume low temperatures, where 
scattering due to charged impurities is important and 
may dominate, and where electron-electron scattering 
plays a secondary role. In the regime of weak momentum 
scattering studied here quantum interference effects are 
also expected to be negligible. 43 We stress once more that 
a conductivity independent of n and rii was measured 
experimentally by taking n — ► 0. The value obtained is 
characteristic of the strong momentum scattering regime 
and is referred to as the minimum conductivity. - At the 
same time, theoretical research on clean samples finds 
an additional conductivity independent of n and rii. It is 
this latter term, rather than the minimum conductivity, 
that is discussed in our work. 

The outline of this article is as follows. In Section II 
we construct the kinetic equation for bilayer graphene 
taking the quantum Liouville equation as our starting 
point. We discuss in detail the role of a general elastic 
impurity potential in this kinetic equation specific to bi- 
layers. In Section III we apply this equation to the study 
of transport in unbiased bilayers and identify contribu- 
tions to the conductivity cx n/rii and independent of n 
and rii ■ In Section IV we discuss the role of the gate and 
show that, if treated in a naive manner, it can appear 
to yield an off-diagonal conductivity in the absence of a 
magnetic field. We discuss the implications of this finding 
for transport in graphene bilayers and monolayers. 



II. KINETIC EQUATION 

The formalism parallels that used in Ref. HH and its 
exposition below is correspondingly abbreviated. The 
system is described by a density operator p. Evaluation 
of p in the steady state allows one to calculate expecta- 
tion values such as that of the velocity operator. Very 
generally, p obeys the quantum Liouville equation 

Here H is the band Hamiltonian, H E — eE ■ r the inter- 
action with the external electric field E, and U the im- 
purity potential. We project the Liouville equation onto 
a set of time-independent states of definite wave vector 
{|fes}}, in which (ks\p\k's') = p s k s k , = pkk' and similarly 
Hkk' , H kk , , and Ukk' ■ For bilayers the index s runs from 
1 to 4 as will be shown below. We refer to pkk' as the 
density matrix. Matrix elements of Hkk' = Hk 5kk' are 
diagonal in k but off-diagonal in ss' , and similar for H k . 
Matrix elements of Ukk' are off-diagonal in k. Elastic 
scattering is assumed and the average of terms Ukk'Uk'k 
in the disorder potential over impurity configurations is 
ni\Ukk' \ 2 /V: where V is the crystal volume and Ukk' the 
matrix element of the potential of a single impurity, pkk' 
has a part fk diagonal in fc, and a part off-diagonal in k. 
We will be interested in fk since most operators related 
with steady-state processes are diagonal in k. From the 
Liouville equation an effective equation is derived for fk 
in the first Born approximation, valid for spr/h ^ 1 

^ + ^{HkJk] + j(fk) = Z k , (2) 

where the source term = (eE/h) ■ (dfo/dk), the equi- 
librium density matrix fo(Hk) is given by the Fermi- 
Dirac function and the scattering term takes the form 

/>oo 

J(h) = (m/h 2 ) / dt' [U, e-^'/^U, f] e i&t '/% k . 
Jo 

(3) 

The low-energy bilayer graphene Hamiltonian is 3 - 

/ hvke 1 ^ t± \ 

Hvke-^ , , 

k ~ t± hvke~^ ■ [) 

\ hvke^ / 

Along the diagonal are two 2x2 submatrices which rep- 
resent the Hamiltonians of individual layers, in which 
v w 1.1 x 10 6 m" 1 stands for the (constant) Fermi veloc- 
ity of graphene. These layers are coupled by the inter- 
layer tunneling parameter t± « 0.3 eV. Although Eq. (j4|) 
does not include the so-called trigonal warping terms, 
this model captures most of the important physics^ If 
we define = \/t 2 L + 4h 2 v 2 k 2 , the energy eigenval- 
ues can be labeled as Eki — \{^k + t±) = —£k4 and 
£k2 = \{^k — t±) = — £fc3, independent of the direction 
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of k. The energies Eki and Ek2 are positive, whereas £^3 
and £^4 are negative. The two bands Ski and £^3 touch 
at k = 0. We work henceforth in the basis of eigen- 
states of Hk with the eigenvectors labeled \uk s ) so that 



|fes) 



\uf~ s ). In this basis the Hamiltonian is diago- 



nal and has the form Hk — diag(efei, £^2, £&3, £fe4)- Nev- 
ertheless, one must be careful in writing down the kinetic 
equation in this basis. As the basis functions, namely the 
eigenvectors of H^, depend on the magnitude and direc- 
tion of the wave vector k, the ordinary derivative with 
respect to k must be replaced by the covariant derivative. 
The action of this covariant derivative for example on /o 
is given by Df /Dk — dfo/dk — i[TZ, /o]. The connec- 
tion matrix 7t which enters the covariant derivative has 
elements 1Z SS > — (uk s \iduk s ' / dk) . It is easiest to work 
out the derivatives with respect to the magnitude k and 
polar angle 9 of the wave vector, 



IT = 



ihvt± 



/o -1 0\ 

1 

10 

Vo -1 0/ 



/ 1 - 



K = - 



hvk 





V 



lirk 

1 

-1 







-1 
1 

hvk 



\ 



hvk 

1 / 



(5a) 



(5b) 



These expressions will be needed when constructing the 
velocity operator, as well as when determining the source 
term in the kinetic equation in an electric field. Lastly, 
since the Hamiltonian of graphene bilayers and single 
layers does not depend on the spin of the particles, fi- 
nal results will contain a factor of 2 from the sum over 
the spin, as well as an additional factor of 2 which takes 
into account the twofold valley degeneracy of graphene. 
Therefore the final expressions for the conductivity must 
be multiplied by an overall factor of 4. 

The matrix elements of the scattering potential U due 
to a single impurity in the basis of Hk eigenstates are 



(ks\U\k's') = Ukk> Mki = U kk , M ky 



(6) 



where Ukk' are the matrix elements of U between plane 
wave states. The band indices s and s' will be henceforth 
suppressed and quantities such as Mkk' will be treated 
as matrices in the subspace spanned by the four bands 
under consideration (with the band index s the same as 
that introduced above.) The scattering term J(fk) ap- 
pearing in the kinetic equation simplifies considerably if 
we assume that the tunneling parameter t± 3> hvk. This 
assumption is valid for carrier densities up to approxi- 
mately 10 12 cm~ 2 . At these densities only one of the 
bands is occupied: for electron (hole) doping this is the 
band labeled 2 (3). We expand Mkk' in the ratio Hvk/t± 
up to order 1 , and we find that the term of order 1 van- 
ishes identically. We label the incident wave vector by 
fc, the outgoing wave vector by k' and the polar angle of 



the outgoing wave vector k' by 9' . If we define 7 = 9' — 9 
as the relative angle between incident and outgoing wave 
vectors, Mkk 1 has the simple diagonal form 



M kk ' = 



diag(l, COS7, COS7, 1). 



III. TRANSPORT WITHOUT GATES 



The most transparent solution to the kinetic equation 
is found by dividing all matrices in the problem into a 
diagonal part, denoted by the superscript d, and an off- 
diagonal part, denoted by the superscript od. In the case 
of the density matrix the diagonal part ft represents 
the fraction of carriers which are in eigenstates of Hk, 
while fk d is the fraction of carriers which are a contin- 
ually changing mixture of eigenstates of H k - Diagonal 
matrices commute with Hk while off-diagonal matrices 
do not. The kinetic equation is correspondingly divided 
into equations for the diagonal and off-diagonal parts of 
the density matrix, which are coupled by scattering 



r) f d 
f) f od i 



(8a) 



(8b) 



P d and P od are projection operators which single out the 
diagonal and off-diagonal parts of matrices respectively. 
To solve these equations, we search for terms in the den- 
sity matrix of lowest orders in h/{eFT) or equivalently 
those terms of lowest orders in rii. Inspection of Eq. 
(l8al) shows that, due to the absence of the commutator 



[Hk,fk], ft starts at order n~ x while the leading term 
in is independent of rii (in other words order zero) 4^ 
For weak momentum scattering therefore we only need 
to consider the effect of the scattering term acting on the 
diagonal part ft of the density matrix. This reduces to 
the simple form 



p d ju d k ) 



8H 3 v 2 



d6 ' \U 

Z7T 



2 F d (j)(ft-ft'), (9) 



with the diagonal matrix -F d (7) given by 

F d (j) =diag(2,l + cos27,l + cos27,2). (10) 

This matrix is due solely to the overlap of eigenstates 
at different wave vectors. The equation for the diagonal 
part ft of the density matrix in the steady state, in which 
the time derivative can be dropped, reduces to 



J \U kk '\ 2 F d { 1 ) (f d f d ) = St (11) 

The solution to this equation is found simply as ft = 
Y,f,T m , where r m , which plays the role of a momentum 
scattering time, is a diagonal matrix given by 

^ 1 = ^/f^ fc ' |2F£iW(1 - COS7) (12) 
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with the matrix elements r m ii = r TO 44 = t + and r m 22 = 
r TO 33 = t_ . This allows to write the electric- field- induced 
correction to the diagonal part of the density matrix as 



(per valley and spin) 



f d _ 
Jk — 



2ehv 2 E ■ k 



diag(r + Sx,T-5 2 , -t_ 8 3 , -r+ <5 4 ), 

(13) 

where we have used the abbreviation S s — 6(ek s — £f)- 
We would like to stress that, although r m has the units of 
time, its elements do not correspond to actual scattering 
times and their energy-dependence does not come explic- 
itly through any of the band energies, but rather through 
Afc, which represents the difference between band ener- 
gies. The form of Xk implies that the dominant contribu- 
tion to T m comes from t±. Therefore, whereas /j? cx n^ 1 
does indicate that scattering tends to keep the Fermi sur- 
face near equilibrium, an expression for this function in 
terms of a momentum relaxation time cannot be formu- 
lated in bilayers. 

The next step is to evaluate the off-diagonal part fZ , 
Firstly, the effective source term in Eq. (|8b[) contains the 
off-diagonal projection of the scattering term acting on 
/j?, given by P od J{ft)- We find that this projection pro- 
duces corrections of order (hvk/t±) 2 and therefore can be 
omitted. Secondly, when carrying out the time integral 
that determines the correction due to ££ d we allow the 
field to have a time dependence e~ luJt , taking the limit 
lu — ► at the end. (An unphysical negative conductivity 
may be obtained if this procedure is not followed^) /£ d 
is found using the time evolution operator e~ lHt ' fi 



f" 

Jk 



od 



lim / dt' e~ nr ' e ™(t'-t) e -iH h t'/h E od e iH h t'/h 

(14) 

with r\ > a regularization factor. The time integral 
results in a series of 6- functions of the form 5(sk s — £ks'), 
with s =/= s' . The only S- function that can take nonzero 
values is 5{ek2 — £fe3), due to the two bands that touch at 
k = 0. The problem of finding fg d therefore reduces to 
finding its matrix elements in the subspace spanned by 
|ufe2) an d |ufe3). The only nonzero matrix element is 

/fe23 = ~ m& 2k ° }}™{fo( £ k2)-M £ k2-rKv)}s(e k2 -^ 

(15) 

Thus the bands that touch at k = give f£ d , much like 
in single-layer graphene. Finally, the diagonal projection 

pd J{fk d )i which would act 

as an effective source in the 
equation for f£, is also of order (hvk/t±) 2 and is omitted. 

We determine separately the contributions to the elec- 
trical conductivity due to each term in the density matrix 
for electron doping, for which only the band labeled 2 is 
occupied. We require the expectation value of the cur- 
rent operator j = —ev, where v = (\/K)DHk/ Dk is the 
velocity operator in the basis of eigenstates of H^, and 
the covariant derivative D/Dk has been defined above. 
fu, the fraction of carriers in eigenstates of Hf~, yields 



mrH v t 



± 



h (4Tr 2 h 4 v 4 )n 2 + {4irh 2 vH 2 L )n + t 4 



(16) 



with the dimensionless quantity £ = r_Afc/fi.. In bilayers 
the screening wave vector is independent of the number 
density,^* and (an for long-range scatterers and is a 
constant for short-range scatterers. For short-range im- 
purities, at low densities a^ x cx n (Ref. [H) as in single- 
layer graphene but as the density increases nonlinear 
terms in n become more pronounced. For long-range 
impurities a x d x cx n 2 at low densities. The dependence 
of a^ x on t± indicates that a^ x is due to carriers from 
both layers. Nevertheless, the leading term cx tj 2 , im- 
plying that at low densities interlayer tunneling hinders 
the transport of charge. 

For a nonzero chemical potential /i, fg d gives (per val- 
ley and spin) 



ire 
~4h 



lim 



1 



1 



I + e -/3(>/2+W4) 1 4- e -0{fi/2-hu/4) 

(17) 

where (3 — 1/kT. To obtain the dc result at T = 0, one 
must take the limits T — > and u> — > 0; yet the result 
depends on the order in which these limits are taken. If 
T — > first the result is ire 2 / (4/i), whereas if u> — > first 
the result is zero. The same conundrum is present in sin- 
gle layers of graphene.— At present it is not clear whether 
this term is finite, and the theory presented in this pa- 
per up to now does not offer an indisputable solution. 
So far neither theory nor experiment can disambiguate 
this issue. Experiment could provide a conclusive answer 
if clean samples with zero carrier density were available, 
but that remains a daunting task. For weak momen- 
tum scattering, u x od is considerably smaller than a^ x and 
cannot be extrapolated conclusively from a plot of n/n,. 
Fortunately however, the analysis presented in the fol- 
lowing section on biased bilayer graphene shows that, if 
we consider the conductivity due to a gate the answer 
can be found. 



IV. GATE EFFECT ON TRANSPORT 

We have determined so far the steady-state density 
matrix in unbiased graphene bilayers. We study next the 
interesting case of biased bilayer systems, in which the 
gap can be modified by the application of an external 
gate voltage V g . The gate voltage gives rise to an addi- 
tional term H? in the Hamiltonian, which in the basis of 



eigenstates of Hk has the form 



H 9 k = 





( 





—hvk 





tx 


eV g 




-hvk 





-t± 





2A fc 


{ 





-t± 





—hvk 




t± 





—hvk 






(18) 



We treat the gate potential — eV g in first-order pertur- 
bation theory. It is easily checked that to first order in 
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\eV g \ and hvk/t± the contribution of to the scattering 
term J(fk) hi the kinetic equation is zero. The gate gives 
rise to an additional source term in the kinetic equation, 
which takes the form 

(19) 

I 



Several other terms either vanish as lu — > or give zero 
contributions to steady-state expectation values after in- 
tegration over wave vector. These terms are omitted here 
for simplicity and without loss of generality. 

I 



The term <jg V , like and like the term analogous to 
in single-layer graphene,— depends on the order of 
limits. Nevertheless, this term is the only contribution 
to the off-diagonal conductivity, a fact that can clarify 
the correct order of limits. Crystal symmetry and On- 
sager relations imply that an off-diagonal conductivity 
a xy requires time reversal symmetry breaking, for exam- 
ple through the presence of a magnetic field. Therefore it 
must be argued on physical grounds that an off-diagonal 
term in the conductivity such as that found in the cur- 
rent work should not exist 4^ The only resolution to this 
is to take the limit lu — > hrst. These observations imply 
that o~o% in bilayer and single-layer graphene should be 
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Following the same procedure as in the unbiased case 
discussed above, we find that the nonequilibrium correc- 
tion /? to the density matrix due to the gate gives only 
an off-diagonal conductivity (per valley and spin) 
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zero. 

We note that the findings of this section do not con- 
flict with the presence of a minimum conductivity a xx = 
4e 2 /(Trh) in ballistic graphene, which has been measured 
experimentally. 4 - The ballistic regime is qualitatively dif- 
ferent in that the mean free path greatly exceeds the 
sample size. 
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